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Abstract 

This paper is the second part of a threefold article, aimed at solving numerically the Poisson 
problem in three-dimensional prismatic or axisymmetric domains. In the first part of this 
scries, the Fourier Singular Complement Method was introduced and analysed, in prismatic 
domains. In this second part, the FSCM is studied in axisymmetric domains with conical 
vertices, whereas, in the third part, implementation issues, numerical tests and comparisons 
with other methods are carried out. The method is based on a Fourier expansion in the direction 
parallel to the reentrant edges of the domain, and on an improved variant of the Singular 
Complement Method in the 2D section perpendicular to those edges. Neither refinements 
near the reentrant edges or vertices of the domain, nor cut-off functions are required in the 
computations to achieve an optimal convergence order in terms of the mesh size and the number 
of Fourier modes used. 

Date of this version : July 7, 2005 

1 Introduction 

The Singular Complement Method (SCM) was originally introduced by Assous et al. 
[01 Ej, for the 2D static or instationary Maxwell equations without charges. It was then 
extended [U03 to the fully axisymmetric case, i.e. axisymmetric domains and data, with 
or without charges. The SCM has been extended in ^T] to the 2D Poisson problem. 
As noted in ^21) further extensions to the 2D heat or wave equations, or to similar 
problems with piecewise constant coefficients, can be obtained easily. Methodologically 
speaking, the SCM consists in adding some singular test functions to the usual Pi 
Lagrange FEM so that it recovers the optimal H 1 -convergence rate, even in non-convex 
domains. In the fully axisymmetric case, one may simply add one singular test function 
per reentrant edge, and one per conical vertex of sufficiently large aperture. 
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There exist a couple of numerical methods in the literature for accurately solv- 
ing 2D Poisson problems in non-convex domains. The SCM is clearly different from 
(anisotropic) mesh refinement techniques |2L)l I16j . and can be applied efficiently to in- 
stationary problems (see Remark 4.1 of |12|). since it does not need the refinements 
of the mesh and thus large time steps may be allowed. However the anisotropic mesh 
refinement methods have one advantage: they require only a partial knowledge of the 
most singular part of the solution. 

The numerical solution of 3D singular Poisson problems is quite different from the 
2D case, and much more difficult. This is a relatively new field of research: most 
approaches rely on anisotropic mesh refinement, see for instance |16[ I17j and Refs. 
therein. To our knowledge, this series of papers is the first attempt to generalize the 
SCM for three-dimensional singular Poisson problems. 

The rest of the paper is organised as follows. In the next Section, we define the ge- 
ometry of the axisymmetric domain and the suitable framework for the study of the 
Poisson problem in f2 using a Fourier expansion with respect to the rotational angle 9, 
namely, weighted Sobolev spaces over the meridian section uj. This suggests a frame- 
work for building the Fourier Singular Complement Method (FSCM) for accurately 
solving the Poisson problem, using a Fourier expansion in 9, and an improved variant 
of the Singular Complement Method in uj. In Section 3, we study theoretically 
this variant of the SCM, based on a regular-singular splitting of the solution u k to the 
2D problem (JTHEJ- The main feature of the splitting is that it is chosen independently 
of the Fourier index k as soon as \k\ > 2; this independence is important, and very 
helpful, from the computational point of view. Section 4 presents a few results of finite 
element theory in the weighted Sobolev spaces. In Section H3 the SCM is considered 
from a numerical point of view, to approximate u k accurately, via the discretization 
of the splitting. In the Section H3 we build the numerical algorithms which define the 
FSCM, and we show that it has the optimal convergence of order 0(h + -/V -1 ), where 
h is the 2D mesh size and N is the number of Fourier modes used. 



2 Poisson problem in axisymmetric domains 
2.1 Geometric setting and notations 

In this article, we consider an axisymmetric domain £1, ggenerated by the rotation of a 
polygon uj around one of its sides, denoted j a . The boundary of uj is hence Ouj = 7 a U7&, 
where 7& generates the boundary T of 0. Thus, ft can be described as: 

n = uxS l U 7 a . (1) 

The natural cylindrical coordinates will be denoted by (r, 9, z). The geometrical singu- 
larities that may occur on T are circular edges and conical vertices, which correspond 
to off-axis corners of ^ and to its extremities. Figure Q precises the various notations 
associated to these singularities; a more complete description of the geometry of u can 
be found in Hj. 
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Figure 1: Notations for the geometrical singularities; e: reentrant edge; c: conical vertex. 



The problem under consideration is once more the homogeneous Dirichlet problem: 

o 

Find u £ i? 1 (f2) such that 

-An = / in Q, u = on T, (2) 

with / 6 L 2 (Q). Non-homogeneous Dirichlet boundary conditions, or (non-)homogen- 
eous Neumann boundary conditions can be handled in exactly the same manner. 

As will appear in the sequel, the problem (|2j) will be singular, i.e. its solution will 
generically not be in H 2 (Q) — as it would be the case in a regular or convex domain — 
iff there are reentrant edges or sharp vertices in T. Sharp vertices are defined by the 
condition (see Figure ^): 

v c < -i, where: v c = f min |^ > : P v ^cos ^\ = j , (3) 

and P v denotes the Legendre function. This is satisfied iff it/ (3 > tr/ft* ~ 130° 48'. 
From now on, we shall assume that there is exactly one reentrant edge e (of aperture 
7r/a, with 1/2 < a < 1) and one sharp vertex c, and we shall omit the superscript c 
in v c . 



Other notations. We denote by r max the supremum of the coordinate r on u, and 
by ao and a\ two fixed numbers such that 

1/2 < «o < a and 1/2 < ct\ < min(a, v + 1/2). 

We also introduce 2D neighbourhoods tu e and uj c of e and c respectively. They stay 
away from all sides of duo except the two ones that meet at the relevant corner. To 
them we associate cutoff functions denoted n(p), which vanish outside u> e or uj c and 
depend only on the distance to the corner. 
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2.2 Fourier expansions 

The functions defined on f2 will be characterised through their Fourier series in 9, the 
coefficients of which are functions defined on to, viz. 



f(r,e,z) = -= /V,*)^ 

V^vr , 



_ifc0 

/ V > *} > 

k=— oo 

and the truncated Fourier expansion of / at order N is: 



1 N 

fW{r,9,z) = -= V / fc (r,^)e 



fc=— JV 

The regularity of the function / in the scale H S (Q) can be characterised by that of the 
(f k )kez in certain spaces of functions defined over uo jlUl §§II. 1 to II. 3], namely: 

/ G H S (Q), s > <=> Vk G Z, / fc G # ( a fc )M> 

where the iT? fc s (a>) are defined in turn with the help of two different types of weighted 
Sobolev spaces. We shall now give these definitions for the values of s and k chiefly 
needed in this article. The notations for the various spaces are the same as in jlUj . 
where the interested reader can find the proofs and the most general versions of the 
subsequent statements. 

First, for any r G M we consider the weighted Lebesgue space 

Lr(to) = f I w measurable on to : / / \w(r, z)\ 2 r T drdz < oo 



This space, as well as all the spaces introduced in this article, is a Hermitian space of 
functions with complex values. The scale (H^.(u)) s>0 is the canonical Sobolev scale 
built upon L^(lo), defined for s G N as: 

H s t {uj) = \w G L 2 t {uj) : d e r d™w G L 2 T (u), W, m s.t. < £ + m < s} , 

and by interpolation for s ^ N. We denote by || • || S)T and | • | SiT the canonical norm 
and semi- norm of H^{u). 

A prominent role will be played by L 2 (u>); its scalar product is denoted (-|-), without 
any subscript. Upon this space, we build another, dimensionally homogeneous Sobolev 
scale (Vl(u)) s>0 , defined as: 

V?(u>) d = {w G H{(u) : r e+m - s df.d™w G L\(u), \f£,m s.t. < £ + m < [s\} , 

where [s\ denotes the integral part of s. One can check that the general definition 
reduces to 

V{(u) = G H\{u) : d}w\ la =0, for < j < s - l} , 
when s is not an integer; while for the first values of s G N, we have: 

y?(w) = £?(w), v\{u) = h\(u) nL 2 ^), vI(uj) = h\{u) nHl^u). 

The canonical norm of Vf (w) is denoted by ||| • ||| Si i; it is equivalent to | • | Sj i except 
for s G N*. 
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We are now ready to define the most useful spaces of Fourier coefficients. 
Lemma 2.1 The spaces Hf^Ju), for s = 0, 1, 2, are characterised as follows. 

H? k) (u) = L?(w), Vfc, = ffjH, ^jH = V\(u), V|fc| > 1; 

i? ( 2 0) (u;) = {u; G : 3 r ™ G L^H} , flf ±1) M = G flf (w) : w\ la = 0} , 

H? k) (o J ) = V 2 (u J ), V|fc|>2. 

T/ie definition for the other values of s will be given when needed. 

Remark 2.1 T/ie scales H\(oS), V\{oj), and H^Joj) (for any k) can be extended to 
negative values of the exponent s, by the usual duality procedure with respect to the 
pivot space, which is L\(uj) in all cases. 

o 

In order to handle the Dirichlet condition, we introduce the subspaces H\(uj), 

o o 

Vj;(u;), H_i{u) of functions which vanish on 7&. The difference in the notation is 
to remind that the functions of V];(u;) and H_ 1 (u) C V^(a;) automatically vanish 
on 7 a in a weak sense |18| Prop. 4.1]. This difference is of course important when it 
comes to discretisation by Pi finite elements. 

Similarly to the prismatic case, we introduce the anisotropic Sobolev spaces 
h\n) = H 1 (S\L 2 1 (u J ))={feL 2 (n):d e feL 2 (n)}; 
h 2 (n) ^ H 2 (S\L 2 (u;))={feh 1 (n):d 2 feL 2 (n)}; 
they are identical to the H°> s (n) of fUR Eq. (II.4.16)], for s = 1, 2. 

The next Lemma summarises the completeness results whose proofs can be found 
in pm Chapter II] or [Hj. 

Lemma 2.2 The following characterisations hold: 

f£L 2 (n) ^ Vk eZ, f k G L\(u), and: ^ \\f k < oo ; (4) 



k=— oo 



+oo 2 



feh s (n) Wk 6Z, f k G L 2 (uj), and: ^ k 2s f k Q ^<oo; (5) 



k=— oo 

and the canonical norms of L 2 (Q) and h s (Q) are equal to the square roots of these 

o o o 

sums. Moreover, defining = HhJu) C\H\{oS), viz. H\{u) for k = and V\{uj) for 
the other cases, we have: 

o + °° II 2 

/ e H\Q) <=> \/k G Z, f k G V (k) and |/|^ 1(n) = ]T / fc < oo, (6) 

(re) 



fc=— oo 



where the norm \\w\\? k \ = \w\ 2 ^ + k 2 ||u>|lo,-i- 

As we did in the prismatic framework, we define the relation operators < and ~ 
as follows, a < b means a < C b, where C is a constant which depends only on the 
geometry of the domain u, and not on the mesh size h, the Fourier order k, or the 
data / of the Poisson problem, a ~ b denotes the conjunction of a < b and b < a. 
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2.3 Singular Poisson problem in 2D 

Denoting by u k and f k the Fourier coefficients of u and / in (J2J, we see |101 §11.4] that 
for any k, u k is solution to the following singular Poisson problem in u>: 
Find u k such that 

a k def A k *k ■ i a def 2 Id d 2 k 2 

A k u k = -A k u k = f k m u, where: A k = —» + - — + — - -j, 7 

or z r or oz z r z 

u k = on 76 . (8) 

A special role will be played by Aq, whose values are the traces in a meridian half-plane 
of the Laplacian of axisymmetric functions. We remark that the operators have 
real coefficients, hence the real and imaginary parts of the solution to (J7HHJ) correspond 
to the real and imaginary parts of the data. So, in practice, it will be sufficient to 
consider problems with real data and solutions. 

The variational space associated to ((31EJ) is the Vrw defined in Lemma 12.21 The 
variational formulation reads [1U1 §11. 4. a]: 

k \ I tk 



ak{v*,v) = {r\v), VveV {k) , (9) 
where a k is now the sesquilinear form defined by the norm || • \\( k ), viz. 



a k (u,v) 



„ ^_ k 2 _ 

r vii • VD H uv 

r 



duj. 



(In this text, V will always denote the 2D gradient in the (r, z) plane.) 
Like in the prismatic case, we have the following results. 

Lemma 2.3 Let f £ L 2 (£l), and u be the solution to (0). Then (u^) K converges to 
u in H 1 ^), and (Au^l) „ converges to —f in L 2 (Q). 

Proof: Similar to |12[ Corollary 3.1]. o 

Lemma 2.4 Let f G L 2 (Q), and u be the solution to Then dgu E H 1 ^). 

Proof: One may follow the lines of |12l Corollary 3.2], using the a priori estimates 
of Pa Thm 4.2] to check that (d r d e u, r" 1 d e d e u, d z d e u) G L 2 (0) 3 . o 

Besides the variational space, we shall consider, for each Fourier mode k: 

• the natural space, which is the one to which u k belongs, i.e. the domain of the 
operator A k : 

D{A k ) d ^ e V (k) = H\ k) {u) n H\(u) : A k w G ; (10) 

• the regularised space, i.e. the one to which the solution u k would belong if the 

o 

domain Q were regular or convex, namely H 2 k ^(uj) n H\(uj). 
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In |1(J1 Thm II. 3.1], it is established that the regularised space no longer depends on k 
as soon as \k\ > 2; in Theorem 13.21 we will show that the same occurs for the natural 
space. This suggests that the mode 2 can serve as the "fundamental mode" for the 
high-|/c| modes, just like the mode does in the prismatic case. In contradistinction to 
the latter, the modes and ±1 have to be treated separately, with their own singular 
functions. 



3 Regular- singular decompositions in the 2D do- 
main uj: theoretical study 

We now establish the regular-singular decompositions, for the various Fourier modes k, 
of the solution u k to (jJHSJ), which will be effectively used in the numerical method. This 
parallels the work exposed in the companion paper |121 §4]. 

We shall need the following integration by parts formulae. 

Theorem 3.1 For any u, v £ H\{uS) such that Aqv £ L\(uj), there holds: 



f f f dv 

// {u Aqv + Vu • Vv} r duj = / u — rd'y. 

JJw ' Jj b ° V 



For any w £ H^^ui) such that Aqw £ L\{ui), there holds: 

w 



Aqw > rduj 



|VHIo,-i 



2 IM 



0,-3 • 
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(12) 



Proof: Eq. (jllj) is the expression, in a meridian half-plane, of the usual Green formula 
applied to axisymmetric functions. To prove ()12|) . we first note that there holds, in the 
sense of distributions in u: 

w\ \Vw\ 2 2w dw 

But w £ H^^to) implies w £ L 2 _^{uj) [3 Lemma 4.9], i.e. r~ 2 w £ L\(u); so the above 
function is integrable with respect to the measure rdui, and we can apply (fTTj) with 

u = r~ 



Vw ■ V 



2 w and v = w: 

T del 
'1 = 



def 



r 

\Vw\ 2 



— Aqw > r duj 
r du — 2 



r dcu 



1 dw\ l w . 
- -5— I — I r duo 

r or J 



I Viol 



0,-1 



21*. 



Now, we treat I2 by the usual integration by parts formula of order one: 



dw 
dr 



w 



+ w 



d_ 

dr 



w 



duj 



dw 



w 



w 



dr r 



dw 
dr 



+ u>r — o 



W A 

w — v r doj 

2 

3 



r 



du) 



2K/ 2 -2||to| 



0,-3 






0. 
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Hence, K/i = ||Vio||o _i — 2 ||io||q _ 3 . o 
3.1 Modes \k\ > 2. 

From jlUl §11.4], we know the following facts. The solution u k to (JTHHJ) is regular 
everywhere except in the neighbourhood of the reentrant edge, and it can be written 
as: 

«* = «*+ a* ss, with: k€flf fc) ( w )n^( w ) = v?( w )n^H > (13) 

\S%M) = r,(p)e-\ k \<> p a sin(a0). 
As a first consequence, we have the following 

Theorem 3.2 Let w £ D(A k ). Then: 

o 

• w has a V\ regularity near the axis, hence w £ H i _ 1 (uj) C L 2 _^(uj), and both Aqw 
and r~ 2 w are in L\{uj). 

• w has an H 1+a ° regularity near the reentrant edge, so its global regularity is w 6 
V\ +a °{ijj); and there holds: ||iu|i_|_ a0! i < C(k) ||Afc«;||o,i. 

In close analogy to the orthogonal decomposition of L 2 {ui) introduced by Gris- 
vard T4J p. 45], we have: 

L\{u) = A 2 [H 2 2) {uj) n V\(lo)} N 2 , (14) 

where N 2 is a space of singular harmonic functions defined by 

N 2 = {p € L\(uj) : A 2 p = in oj, p = on each side of 7^} . 

Here, as well as in the subsequent definitions of N\ and No, the boundary condition 
on the sides of 7^ is understood in the suitable space, which is the trace in a meridian 
half-plane of the space H(Ti) defined in [2>, Definition 5.4]. Following the same line 
of proof as in [3J §3], it is not difficult to establish that the dimension of N 2 is equal 
to the number of off-axis re-entrant corners in u, i.e. in our case dim ^2 = 1, and 
N 2 =sp&n{p 2 }, where p 2 can be chosen as: 



2 =S e +p 2 with: \pI^V\(u), _ (15) 
I S e (p, (/)) = Tj(p) p a sm(a^). 



Similarly to §4], we define ip 2 as the element of V];(a;) which solves the Poisson 
problem 

-A 2 (Ps=p 2 s in u. (16) 
Then by the decomposition (|14jl. we can split the solution u k to (0JHJ) as 



u k = u k + c k <p 2 s , (17) 
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where u k £ H? 2 Jlj) n V\(u) = H%duj) n H\(u), and is called the regular part of 

u k . How is this decomposition related to (fT3)) ? Applying ((T3|) to ^ itself gives: 
(Ps = (p 2 + <5 2 5|; observing that all the S| have the same principal part, we deduce 



— u k + Ck<p 2 R , and = Ck 5 2 . Then, using the orthogonality relation (pH)) we infer: 



u R „ , 



II 2 II 2 

^2 Ik»IIo,1 



(-A 2 Sf|p 2 )' 

Calculating this scalar product is rather tedious but can be done using (|15j) and — 
modified so as to avoid the singularity. We find: 

— = 5=— \\p a \\ 01 , (18 
Cfc an u ' 1 

where o = r(e) is the distance from the reentrant edge to the axis (see Fig. 

The following lemma summarises some a priori estimates on u k and q,. 

Lemma 3.1 Let u k be the solution to the Poisson problem Q-EP, then we have the 
following a priori estimates: 

fc 2 ||u fc ||o,-i<r max ||/ fc || ,i, fc|u fc |i,i < ^ll/'llo,!, (19) 

{k 2 -2) ||u fc || ,-3<||/ fc ||o,i, [k 2 -2) 1/2 \u\_ 1 < -^ll/ fe ||o,i, (20) 

||A n fc ||o,i < 2||/ fc || ,i, (21) 
M < k a ~ l ||/ fe || 0> i (22) 
MIlW « lAw S ^^ll/llo,!. (23) 

Proof: The variational formulation @ with v = u k gives: 

l« fc |?,i + ^ 2 IK fc Ho,-i < ll/ fc Ho,i \W k \\o,i < r max ||/ fc || ,i ||u fe ||o,-i , 

this proves the first estimate in (|19|) . Then applying the Young inequality, we further 
obtain 



I k\2 i ; 2 II kn2 max 
|U |ii + fe \\U \\o -i < —T— 



k' 2 



\\f k \\L + ^\\A\U 



which leads to the H\ semi-norm estimate in 1)19(1 . Similarly, multiplying (JJJ) by 
and using ((T2() yields: 



2 7rrfc 



r W 



\u k \l_ 1 + {k 2 -2) Wu%- 3 <||/ fc ||o,i||n fc || ,-3, 

and we obtain the two estimates in l(2U(l by a similar reasoning. Then ()21() immediately 
follows from Aou k = f k — k 2 r~ 2 u k . 

The formula 1)18(1 implies: |cfc| ~ |Afc|; thus, the estimate (|22j) is clearly equivalent 
t° : |Afe| < This, in turn, can be obtained by following the lines of [141 §2.5.2] 

or [1J §5.1]. As a matter of fact, the latter reference shows that, away from the axis, 
the weights in the Sobolev spaces and the exact form of the modified Laplacian under 
consideration are of no importance. 
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Now, setting /* d = -A k (u k - c k tpf), i.e. 

fn = - A on k R + ^n k R , 



k 2 



(24) 



one concludes, like in the above references, that ||/^|| 01 < ||/ fe ||or Expanding the 
squared norm of the equality (|24|) and using ()12l) then yields: 



A < 



0,1 



+ 2r 



u , 



1,-1 



+ (A: 4 - 4fc 2 ) 



0,-3 



/ 



0,1 



(25) 



On the other hand, there holds: u R G i?^(cj) H iif];(u/), and within this space the 
canonical norm of H 2 Q ^(u>) is equivalent to the norm ||Aou>|| i @J Lemma 4.7]. So, we 
have both \u R \ 2 1 < ||/ fc || j an d J x 1 < A;^ 1 ||/ fc || and we obtain by interpolating 
in the scale H\(u) that: |w^| 1+qo l < k a °~ l ||/ fc || r We then derive by adapting 
the proof of Lemma 4.1 of ^2]- ° 

Lemma 3.2 The regular part u k and the singularity coefficient c k in are given 
as the unique solution of the coupled system: 



a k (u k ,v) + c k a k (p 2 s ,v) = (/ fc | vj , G V\(u) : 



\\p 2 s fo,i + » [|v2ll,-i + 2||^llo,-3]) C fc + Ai 



U k \Ps 



(26) 
(27) 



where the symbol fi = f k 2 — 4. And u k and c k have the following stability estimates: 



l^ll(fc)<fc||/ fe Ho 



Icfcl < 2 



ll/ fc ||o,i 

l|Pal|o,l 



2,1 



< 



10,1 • 



We omit the details of the proof, which is very similar to that of Lemma 4.2 of |12j . It 



makes use of the result: 



- I^sli -i + 2 ll^sllo -3' w hich directly follows 
from H12|) and (|16j) , as p 2 and y 2 , are real. The representation formula for the singularity 
coefficient is: 



Ck 



{f k -{n/r 2 )A^f k 



|| p 2,|2 



(28) 



s 110,1 



The scalar product (r 2 A k 1 f k \ p 2 ) in (|28[) is defined thanks to Theorem l3.2l We shall 
see — and this will be of practical relevance — that it can be written as (A k ~ 1 f k \ r~ 2 p 2 ) , 
i.e. p 2 G L 2 _ 3 (oj). This is a consequence of the following lemma. 

Lemma 3.3 The dual singularity p 2 is of V 2 regularity near the axis. It admits the 
following splitting: 

P 2 s=P 2 +P 2 p , P 2 £V\(lu), p 2 p = {£j 2 P~ a sin(a0). (29) 
Similarly, the primal singularity (p 2 can be represented as: 

2 



£ 2 + 5 2 <p 2 p , (p 2 G V\{u), V>1 = (£) P a sm(a ( j ) ). 



(30) 
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Proof: Let < a" < a' < a; we consider a cut-off function x such that x( r ) = 1 f° r 
r < a" and x{ r ) = f° r r — a 'i as wen as the domain u' = {x G : r(a;) < a'}. This 
domain has no off-axis reentrant corner (see Figure ^) , so there are no singularities 
of A2, either primal or dual, in u>'. 

As we stay away from the reentrant corner, the splitting 1)151) shows that p 2 G 
l^}(o/). Thus, xPs ^ ^lO^O an d it vanishes on du>' . Moreover: 

A 2 (x rf) = X &2P 2 S + V X • Vp 2 + p 2 . A X + & ^ G L? (a/) , 

since the first term is identically zero, and the other three are smooth and vanish near 
the axis. We conclude from Theorem 13.21 and the absence of primal singularities, that 
XP 2 S G V\(uj'), i.e. p 2 is V 2 where x = 1. 



Now, using (see Figure P): 

T\ 2 (r — a){r + a) 2 



a J a? a 



p cos(0 + 4>o) + h.o.t., 



we remark 

S e - v\ = P l ~ a (si(<£) + h.o.t.) G H\u e \ SI -v 2 p = p 1+a (g 2 ((f>) + h.o.t.) G H 2 (u> e ), 

since the functions gi^i^ 1 ) as wen as the higher-order terms (h.o.t.) are smooth. More- 
over, thanks to the factor (r/a) 2 , p 2 p and ip 2 p are of V 2 regularity near the axis. The 
smoothness of these functions in the rest of the domain yields S e — p 2 p G ^(cj), 
S 2 e -^ G V\(uj). This proves and (3UJ). o 



3.2 Modes k = ±1. 

o 

As we can see from Lemma l2.ll the variational space is still y}(o;); but the regularised 
space has changed. Once again, the only singularities are located at the reentrant 
edges. Hence, the solution u k to (JZHHlli with k = ±1, can be split as: 

u k = u k R + x k s^ with: (<eFf ±1) Hnf{H = FfHnf{H, (31) 

yS e ±1 {p,(t)) =r/(p)e p p a sm(a(f>). 

As a consequence of Theorem 13.21 ip 2 G D{A\), and the decomposition 1)17)1 is still 
valid in this case. However, that singular function belongs to a space which appears 
too constrained for the modes ±1: it is even better decaying near the axis than the 

o 

functions of H\{uj)f\V\{u); moreover, this decay is lost in the discretisation by Pi finite 
elements. So the representation formula ()28|) . though valid at the continuous level with 
p = —3, is numerically hardly stable and its use would deteriorate the convergence rate 
of the SCM. 
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So, it is better to use singular functions that are adapted for these modes. Let pi 
be a basis of the dual singular space 

JVi = {p £ L\(uj) : Aip = in ui, p = on each side of 75} , 

and <pl = A^ l p\. These functions were defined and studied in jSJ §4.1], and a numerical 
method was defined. We will introduce below ( fl5.4j) a slight modification of that 
method in order to improve the convergence rate. For the moment, we recall that the 
function u k admits the splitting 



u k + c k (fl, 



(32) 



where u k G HlJu) n H\(oj) = H\{u) n V\(u). As we are in the "usual" SCM 
framework jllj . we have the simple representation formula 



llri llii 



and the regular part satisfies: 

at (u k , v) + c k a x {if] , v) = (^f k \ u) VveV\(u), 
From the above considerations easily follow the estimates: 



u k 


< 






1,1 ~ 




u k 




< 




l+ao,l 





0,1 



0,-1 
\U 



< 



< 

1,1 ~ 



A u k 



< 



0,1 



0,1 



0,1 



\Ck 



< 



0.1 



(33) 

(34) 

(35) 
(36) 



3.3 Mode k = 0. 

© 

Now, the variational space is V( ) = H\(uj), and the regularised space is H? \(lj) n V"( ), 
with H'Lx (lo) given by Lemma 12.11 Moreover, there is one singularity per reentrant 
edge and one per sharp vertex, see [TU1 §11.4] or H §4-4]. The splitting of u° with 
respect to regularity thus becomes: 



u° R +\ e S% + \ c SS, with: < 



u R £Hf 0) (u)nH\(Lu), 

s o(p,4>) = V(p)p a sm(acf) 
{S c v(pA)=ri{p)p u Pv{co S , 



(37) 



Once more, there holds: H? Ju) C H"LJu); yet, once more, attempting to represent the 
singularity at the edge with the help of the function (p 2 s would imperil the convergence 
rate of the numerical method. As a consequence, we shall use the "natural" singular 
functions for this mode (see |3J[!)] for details). The dual singular space is 



Nn 



{p G L\(uj) : Aqp = in u, p = on each side of 7&} ; 
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it is of dimension two, with the basis |ps' e , Ps' c | given by 



P t = T +P r) with: 1 (38) 

__ . o 

and S e is as in ()15|) . The primal singular functions (p s ' J £ -f^i(^) (i = e, c) are of course 

defined as: (/?s'"' = ^o^Ps 1,7 - 

Lemma 3.4 The primal singular function (p^ admits the splitting 

pOj = ^o,j +6 o,j s ^ where: p° R j £ Hf Q) (uj)r)H\(u}), and: (39) 

6 °' C = \\Ps C \\a,i {( l + 2u ) J Pu(cos(f>) 2 sin0d0l , (40) 

S °' e = Vs'Xa/^)- ( 41 ) 
Equivalently, p^ and tp^ can be represented as 

^ eH ?M, (43) 

Proof: Let us examine first the conical singularity: j = c. The splitting (|37|) applied 
to 99s' c yields: </?s' c = y?°' c + 5®'° Sq + 5c ' c Sq. Using the definitions of pi' c and ip°s C we 
deduce 

llrf'lo,! = # c (-Ao^o e I rf' c ) + # e ("Ao5 c | p°' c ) . (44) 
Rewriting the first scalar product in (|44[) as: 

(-A S e | pY) = (-A S e | p° R c ) + (-A S e | S c ) , 

we notice that the second term is zero by a support argument. To evaluate the first 
term, we remark that —Aqp^ c = AqS° vanishing near the conical point, is smooth 

o 

and belongs to L\{<jj). So, both Sq and p°^ c are functions in H\(uj) with Laplacian Ao 
in L\{uj)] and we use twice to obtain: 

(-A 5 e | p /) = (S e | -A p /) = (5 e | A S C ) , 

which again vanishes by a support argument. Finally, the last scalar product in (|44jl 
can be computed e.g. as in ^H] to obtain 



{-A S c \pY) = {l + 2v) / P t 

Jo 



7T//3 

u \COS (t>) 2 sin 



and (|39H4(J[) are proven. This immediately yields (|42H43j) . Similar computations are 
carried out for the edge singularity; 8°' e is computed like 5 2 above. o 
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Lemma 3.5 The solution to —AqvP = f° can be represented as 



u 



- u° + c , e tp° s ' e + c , c f° s ' c , where: u° G H? 0) (u) n ^l(w). (45) 



(0) v 



The CQj are given by the representation formulae: 



(/° I P°s' j ) 



Wps llo,i 

Proof. As the space of singularities is of dimension two, it is enough to exhibit two 
linearly independent functions to have a basis. This is obviously the case of <pi ,e 
and <^s' c , which proves (|45jl. Taking the Laplacian — Ao of this equality and the scalar 
product by pi' c yields, thanks to the orthogonality property: 

(/° I P°s' C ) = C0,e ("Ao^' e | PY) + C0,c \\p° s ' C \\h- 

Then, using the decomposition (J3HJ), we obtain: (^—A (p®' e | p®' c ^ = 5°' e (^—AqSq \ p®'°^J, 

which is zero as seen in the proof of Lemma 13.41 Hence (j46fl for j = c; the case j = e 
is treated similarly. o 

Let us state without proof the elliptic equation satisfied by u°: 

a (u°, v) + c , e a {<p° s ' e , v) + c , c a (^' c , v) = (/° | v) , Vf G H\(u), (47) 
and the stability estimates on the various terms in (|45|) : 

\Ai,i^\\f\^ \Ai +ai ,^\\f\^ ( 4 8) 
l^ |i.i^l|/ Hoi' ko,j] < ||/°||oa- < 49 ) 



4 Interpolation and projection operators 

We consider a regular triangulation of the domain u>, with mesh size n. The space 
spanned by Pi finite elements on this triangulation is denoted V h ; the subspace of 

o o 

functions which vanish on the whole of duj is Vq = V h C\V\(uo)] while = V h nH\(u) 
is the subspace of functions which vanish only on 75. We introduce the usual Lagrange 
interpolation operator Uh as well as the weighted Clement operator Ph. The latter — 
identical to the operator called 11° in §4] — is a local projection operator onto Pi in 
the L\ sense, which does not take into account the nodes of the triangulation which 

o 

stand on dto. Hence, it maps ^}(w) onto V^. 

We now prove a few results on these operators, in the framework of weighted Sobolev 
spaces of fractional order. We begin by a useful density lemma. 

Lemma 4.1 V\ (w) n V\ (u) is dense within V\ +a °{u) n V\{u) = H\ +a °(uj) n V\(u). 

O 

Proof: Let w G V\ +a °(uj) n V\(lj) and e > 0. The construction of w G Vf(u;) n 

o 

V}(u;) such that — u3||i+a ,i < e is decomposed into two steps. The first one will 
not be effectively used in this article, since we do not consider arbitrary functions in 

V\ +aa {uj) n V\(u), but only those which belong to D{A 2 ). 
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STEP 1: From H\ +a °{uj)r\V\{uj) to D(A 2 ). From [TQ1 Thm II.4.8], we know that 

o 

is an isomorphism from Hfj~^ (u>) to Hfj^ (uj) n H\(u) for s < a. Let C(s,k) be 
the norm of this operator. Setting g = — /S. 2 w G iJS~ (w) = H^ ~ (uj), we use the 
density of Hi(u) within H"°~ (u) to construct g' G L\{uj) such that \\g — </||i_ a0i i < 
e/(2C(o!o,2)). Then u/ d = A^V G D(A 2 ) and satisfies ||w - w'\\ 1+ao ,i < e/2. 

O O 

STEP 2: From D(A 2 ) to V\ (w)flFj(w). There remains to find w G ^(u^nF^) 

o 

such that — u>||i +a0j i < e/2. Since -D(^2) = ^i(w) n V^(u;) © span5|, this is 

obviously equivalent to find S G ^f(w) n Vi(cj) such that ||5| — S||i+a 0) i < e' , for 
arbitrary e' . 

We claim that S = S 2 — St does the job for large enough. As a matter of fact, 

o 

S 2 - Sl = p 1+a sin(a0) + h.o.t. G H 2 (uj e ) n H l {u) e ); thanks to the cutoff function rj, 

o 

this gives S 2 — 5| G ^i(w) H Then, it is enough to check that 

S e k -» strongly in F 1+a °(cj e ) as fc -» +00. 

This is done like in |14| Lemma 5.3.3], by using the Sobolev imbedding H 1+a °(uj e ) C 
W 2 ' p (uj e ) with p = 2/(2 - a ); indeed, one calculates II^H^i+^o^) < HSf H^,^) < 

^.2(ao-a)/(2-ao) < 

o 

Lemma 4.2 For any w G (w) n V\{u), there holds: 

\\w ~ P hMh,i £ h ao |||«;|||i +ao ,i, \\w - Phw\\ 0>1 < h 1+ao |||«;|||i +a0l i. (50) 

o 

Proof: Assume first that w G V\{oj) D V}(u>). From [SJ Thm 4.2], we know 
h~ \\w - Phw\\ ,i + \\w - Ph.w\h,i ^ h 3 ^ 1 IHIki for s = 1, 2 ; 

and, from Prop, l.e.l], that y} +a °(u;) is the interpolate of order «o between V^(cj) 
and V\{uj). Hence the two estimates in (|5U|) by a standard interpolation argument. 

o 

Then one generalises to 10 G y i +a °(u;) Pi V^(a;) by Lemma 14. II o 
Lemma 4.3 For w G H\(uS), there holds: 

\\w - Tl h w\x- a ,x &\w- II/Hi-^i < h 1+ ° \w\ 2 ,i. (51) 

for any a G [0, 1] . 

Proof: It stems from |181 Prop. 6.1] or jSJ Prop. 4.1] that 

h~ l \\w - ILywHo,! + \w - U h w\i,i < h \w\ 2 ,i- 
One concludes by interpolating in the scale H\ (u). o 
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5 Discrete formulation, SCM 

In to 15.31 the superscript 2 in p 2 , ip 2 , 5 2 , etc. will generally be omitted. 



5.1 Approximation of the dual singular function p 2 s . 

We start from the decomposition Q29[) . p is characterised by the three conditions 

peV\(uj), p = -p p on 7, -A 2 p = A 2 p p in uj. 
A direct calculation shows that, denoting (ft' = (ft + (fto (see Figure P): 

" a - 1 sin(a0 + <£')■ 



A 2 p p = — ^ ap 



This function is of H 1 regularity near the reentrant edge, and smooth elsewhere, so it 

o 

belongs to the dual of 1/J(u;). However, it should be noticed that A 2 p p never belongs 
locally to L 2 (u> e ). This phenomenon causes the local regularity of p to be weaker than 
in the prismatic case, and dramatically deteriorates the convergence rate of the SCM. 

This inconvenience can be overcome by enriching the principal part with the next 
term in the expansion of p s near the reentrant corner. To do so, we look for a func- 
tion in yj(u;), with a V\ regularity near the axis, vanishing on duj e D 7, and whose 
Laplacian A2 is locally "almost equal" to A 2 p p . First, we look for a local variational 
solution of 

— AQ = a p~ a ~ x sm(a(f> + in u e , Q = on duj e n 7. 

By separation of variables, we obtain: Q(p, 4>) = \ p l ~ a cos 4>' sin(a^) . As the exponent 
1 — a > 0, Q does belong to i? 1 (w e ). To obtain the V\ regularity near the axis, we 
have to multiply it by (r/a) 2 . A simple calculation shows that: 



An 



2 

-) Q 

a 



+ 



- up 

a 

5 r 



-a-l 



sin(a(f) + <f/) 



2d 1 



[sin( 



Hence, the new decomposition: 



Ps 



def 



Pp+p, where: p p = p 



r\ 



a \aJ 



2 

) Q 



a cos 1 



1 



5p 
2a 



sin (a 



COS ( 



<t> + 4>')] 



p~ a sin(a0), (52) 



enjoys the following properties. First, p = p+ (5r 2 /a 3 )Q £ V\(lu) and it vanishes on 
du e n 7. Then, using 

.2 r r (r — a) r ,, 
a a 2 



-)' 

a J 



p COS ( 



we obtain 



v p = A 2 p p = ^3 p 



25 35 

•— sin(a^)) + — a cos (ft' s'm(a(ft + (ft') 



€ (53) 



As — A 2 p = 'dp, we infer by localisation that p E H 1+a ° (uj e ). Elsewhere, the smoothness 
of -dp implies that of p, so p G V\ +ao (u>), and is V 2 near the axis. 
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Now, we are ready to derive the FE approximation of p s . The variable p solves the 
variational problem: Find p G V}(u;) such that 

^ ^ 1 

p = s on dio, and 02 (p,v) = [ti p \ v ) \/v G V^w). (54) 

Similarly to the prismatic case |12l §5.1], we introduce 

• the boundary function s which is equal to the trace of — p p , hence is zero on the 
two sides that meet at the reentrant corner, and smooth elsewhere; 

• the smooth extension s G H\{oS) of s into u; 

• the variable p° = p — s. 

o 

In the variable p°, the problem (|54"|) reads: Find p° G l^i(w) such that 

o 

«2 (p°,v) = (ti p I v) — a 2 (s,v) Vd G ; (55) 

o 

and we have p° G V 1 °(co) n y}(u). Here, too, we approximate p by p^ 1 = UhS + p° h , 
and p s by Pg = p p + II/jS + p^, where solves the approximate FE problem 

a 2 (f£, « h ) = (#£ I « h ) - a 2 (n h S, u h ) Vv h G F h . (56) 

The notation (i?^ | v^j stands for an approximation by a quadrature formula of the 
integral ff^ ti p (r, z) W^(r, z) r dr dz, with i? p (r, z) given by (|53j). As i? p G -^i(w), we can 
suppose that the error caused by this quadrature is bounded as 

< Cj h qi \\w h \\i,i, \/w h G V h , for some C£ > and q x > 1. (57) 

This can be done e.g. by using a sixth-order Gauss-Hammer formula |21l p. 201], with 
seven points inside each triangle, which does not require the unbounded value of ti p . 
Of course, if Wh vanishes on j a , one can replace Hif/Jl^i with the stronger norm ||w^||(2) 
in l[57|). 

Lemma 5.1 Assume qi > 2; i/ien u>e have: 
\\\Ps-p h s \\\ 1A <h a \ 

Proof: Subtracting 1)56)1 from 1)55)1 yields: 

02 (p° - P° h , v h ) = (ti p - I v h ) - 02 (a - n h s, u fc ) V« ft €K h - (58) 
With Vh = p° h — Wh, this implies: 

||p -^||^ 2) > ||p -^||J 2) + 2a 2 (s-n h s,^-^) + 2 (tf p - tf£ | p° - w,,) . 



(tfj - ti p I w h ) 



- Ps 



n.i 
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Now, we set Wh = PhP° ■ Using l|57|l. we obtain 
\\P -Pl\\\ 2) < \\p°-P h p% ) + 2\\s-U h s\\ i2) (M-p\ 2) + \\p°-P h p 
+ Cj/» qi (\\ P l-p°\\ m + \\ p °-P hP ° 
With the Young inequality, the above estimate becomes: 

lb°-^||? 2) < \\p°-P hP % } + 5 \\~s-n h ~ s \\ 2 {2) + 1 11 ' ; :i|2 ■ " ^ 



+ 



P ||(2) + IIP -^P°II(2) 

(2C^ + 1) + '" ft 2C V' (2) + \\p° - P h p% } 



4 

P H(2) 



Thanks to the equivalence of norms || • ||(2) ~ ||| • |||i,i, we are left with the estimate: 

||p° - pllll, < \\\p° - PhP° Hi, + \\\s - n^l^ + h 2 ^. 

By ^1 Prop. 6.1], there holds: |||s — x 1 < | s 1 2,1 ; by Lemma 14.21 we have 

— -P/iP 1)] x a < h ao ||p°||| 1+Q y s an d P° depend only on the geometry of u so 
all their norms can be seen as constants. Hence, as soon as qi > 1, there holds: 
\\\P° -Pllll^h 200 . Finally: 

11^ - pJlli,! = ||p° + ~s-pI- ^u~s\\i,i < lb - Phlli.1 + III s - £ hao - (59) 

The obtention of the L\ norm estimate also follows the prismatic case closely. Here, 

o 

we define w as the variational solution in V}(u;) to 

o 

a 2 (w, v) = (p° - p° h I v ) , \fveV\(uj). 

By elliptic theory QUI ThmII.4.8] we know w £ H^ a °(cj) = V\ +a °{uj) and |H||i+ Qo ,i < 
||p° — p^||g r Its FE approximation Wh solves 

a 2 (w hl v h ) = (p° -p° h \v h ), VveV h , (60) 
so 1 1 1% 1 1 f 2 ) ^ Hp° — P/lllo l' by usrn S Cea's lemma and Lemma H~2*l we infer: 

II™ " V>h\\(2) £ h ao Ml+ao S ha ° l|P° -PaIIo.1 • 
Then, using successively (foTj|) and (|58j). we obtain 

\\P° -P/illo.i = a 2 (w - w h ,p° -p° h )+($ P - $p I w h S j+a 2 (s-U h s,w- w h )-a 2 (s-U h s,w). 

This is bounded by the Cauchy inequality and (|57|) . as well as the duality argument in 
the scale Vf (w): 



< llp°-py ( 2) Ih-fhll 



1(2) 



/i qi 



(2) 



+ P - n ft s||( 2 ) 1 1 to - w h \\(2) + III s - n ft ,s|||i_ a0j i |||u;|||i +a0i i 
< \\P°-Ph\\o,i {h a »h a » + h^+h\s\ 2 ^xh a °+h l+a ° |S| 2il x^ }, 



where we have made use of |18l Prop. 6.1] and our Lemma 14.31 In order to get the 
h 2a ° estimate, we have to suppose qi > 2. Using once more (|51|). we obtain: 



Ps-P 



(61) 



We are also confronted with the task of approximating q s = p s /r 2 . The scalar 
product [z k | q s ) (see (fTl"]) below) is needed to compute the singularity coefficient. 
However, since p$ is just element-wise linear, it is locally proportional to r in the 

triangles which have one or two vertices on the axis; so q$ = f r~ 2 p h s ^ L 2 (uj). This is 
why we cannot hope to control any such thing as \\q s — q^ || Q r 

Yet, thanks to the bounds (|59j) and 1|61|) for p 1 ^ —p s , we do have the weak estimates: 



resp. 



)| < h 2ao |M|o,- 3 , V»6l J _3(w), 
)| < /i Q0 \\v\\q,-i, VveL 2 ^). 



(62) 
(63) 



5.2 Approximation of the primal singular function (p 2 s . 

We start from (|3U|) . which is sufficient for obtaining error estimates similar to those of 
the prismatic case. Using (fTo^l . we see that <p, satisfying (p = —5ip p on du, solves the 
variational problem: 



where: 



a 2 (<p,v) = (p s \v) + 5 (ip p \v), \/v e V^lo), 
— ^ a p a ~ L sin [(a — 1) 



, def . o / a _i 

ipp = A 2 <p P = —ap a 
a 



(64) 
(65) 



We propose the following finite element approximation of (p in V h : 

<fh = ShKhVp + pi, 

where: TThfP is a simple lifting of the boundary condition, cf. |12| Eq. (40)]; the 
singularity coefficient 5h is computed using 5h = — / / (p^) rduj; and pP h G is 
such that iph is solution to the problem: 

a 2 (<fh, Vh) = (p h s | + Sh {f p | «&) , Vfft S Uj 1 . (66) 

Like above, we assume that the quadrature formula denoted by | v^) satisfies: 

(fi, ~ I w h) < Cq h q2 \\w h \\i tl , Mw h £ V h , for some > and q 2 > 1, 

(67) 

where one can replace ||iy/i||i,i with ||w/ l ||(2) if w h\ la = 0. Then, we propose to compute 
the finite element approximation of ip s as: 



<ptl + <WP ■ 
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Lemma 5.2 The following error estimates hold: 

ll^s - III 1,1 %h, \\<ps -fs\\(k) <kh. 

Proof: We follow the lines of the proof of Lemma 5.2 of the companion paper |12| . 
taking care of the extra error generated by the quadrature. Subtracting Q66JI from (|64j). 
we obtain 

a 2 {(p-(ph,v h ) = (ps-Ps I v h )+(6-5 h ) (ip p | v h ) + 5 h U p - ip h p \ v^j , Vv h e V* 
So, for any Wh G V h satisfying Wh — (fh £ 



11^-^11(2) 



< 



11(2) + 2{ 



0.1 



110,1 



+ \S h \ C 2 Q h* (||£-<^|| (2) + ||£ 



< 2 



^ll^ + ^lb-^ll 2 



+ C 



1(2) 



(2) 



+ |5"^| 2 ||Vp|loi + l^| 2 ^ 2q: 



(68) 



But HV'pllo i ^ s a constant of the domain, and the error on the singularity coefficient is 
bounded as 



\6-6 h \ = — 
an 



Hp. 



s 110,1 



P h s 



Ps-P k s 



< h 2ao 



0.1 



hence 



1. With Lemma IS~T1 (|6*S)) becomes 

lb-^||? 2) <lb-^|| 2 2)+^°+^ 2q2 - 



(69) 



(70) 



To obtain an h 1 estimate, it is thus sufficient to assume q 2 > 1. We then derive 
from lj7U|) that, with Wh = Sh^lh<p/S 



^ll 2 2) < h' ao + \6\~ 2 {\s - 5 h \ 2 w\\ + |j,| 2 \\<p - n 



(71) 



As (p £ H\(lo), we have from [18, Prop. 6.1]: \\(p - II/j^|| (2) « |||^ - IL^I^ < h |^| 2jl , 
which with (|B9*|) gives: |||^3 — <Ph\\\i 1 ~ ll¥> — ¥>h||(2) ~ ^> an d finally: 



llkfl-^lli,! < 
Finally, the estimate on 



ti\ 



$ ~ ^|||i,i + I 5 ~ <^l II^pIIi,! i$ ^ 
follows from 



(fc) 



2 

(k) 



<Ps ~ (ft 



(2) 



0,-1 V 4 



(2) 



< A; 2 /i 2 . 
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5.3 Approximation of u k and Ck in decomposition ( 1171) . for 

\k\ > 2. 

Noting that u k and c k solve the coupled system (|26H27j) , it seems natural to formulate 
their finite element approximations as follows: 
Find u\ G and c\ G R 1 such that: 

a k (u k h , v h ) + c h k a k (<p h s , v h ) = (/ fc | v h ) Vv h G V Q h , (72) 

+ [lrflf,-i + 2[|^||g,_3]) 4 + /, (figlrf) = . (73) 

However, like any function in V^, cfg does not necessarily belong to H^^uj) or L^_ 3 (uj). 
This is no serious problem: like in the prismatic case [121 §5.3], we shall rather discretise 
the representation formula Q28|). which we rewrite as follows: 



1 



Ck 



INIo,i L 



f I p« 



p \ z 



(74) 



where t/ s = p s /r 2 and = A k 1 f k G V^(a;) solves 

a k (z\v) = (f k \v), V«€V-i(a;). 

So, we state the 

SCM Algorithm for finding ufc G V/ 4 and cjel 1 . 
Step 1. Find z\ G F h such that 



a k (z k ,v) = (f k \v) Vv 



eV n h . 



Compute c\ as follows: 



c k 



I II 2 

IPs II 0,1 



/* I P 



M I 4 I 5s 



if k<C*h 2 — o 



for some fixed constant C*, and 



4 = if k>C*h 2 — o . 



Step 2. Find u£ G such that 



a k [u h ,v ) +c k a k [ip s ,v 



f k \v) 



G V r 



h 



(75) 



(76) 

(77) 

(78) 
(79) 



Lemma 5.3 For the solution z k to the problem J75| ) and its piecewise linear finite 
element approximation z k in J 76] ), we have the following error estimates 



z k 


7 k 




z h 


z k 


7 k 




Z h 



7 k 7 k 
Z - Z h 





< 


0,-1 






< 


0,-1 






< 


0,1 





10,1 



h 2a fc 2(ao-l) + ^2 



ll/ fc Ho,i, 



(80) 
(81) 
(82) 
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while for the coefficient c^ in \74\) and its approximation c\ in J77| ), we have 



c k — Cfc 



< 



^2a k 2a 0+h 2 A; 2 ) ||/ fc || ,l 



(83) 



Proof: It follows from (|75|) and (|76|) that 



a fe I - z k h ,z k - z\ \ =a k 



z k ,z k -z k ) = (f k \z k -z k 



This implies 



k k 

z — z h 




z k - z k 


2 

< 


f k 




z k -z k 


^ '"max 






z k -z k 


1,1 




0,-1 




0,1 


0,1 




0,1 



0,-1 



hence (|80(l . Then, using Cea's lemma, Lemma 14.21 Thm 7.1 of j!8j . and the bounds 
(|20j) and (|23|) . we obtain another estimate: 



(fc) 



< 





P h z k 


2 








(fc) 


III**- 


-P h z k 


|| 2 


h 2a 


z k 


2 








l+ao,l 



z k -P h z k 



2 

o,-i 



< 



^2a ^(ao-l) + jfi] fk 



+ (k 2 -l)h 2 
2 

0,1 



2 

1,-1 



Of course, a similar bound holds for any g E L 2 (u;), u; = A^ 1 g and -P/iiu. Thus, the 
estimate (|82|) follows from a duality argument like in |12l Lemma 5.3]. Moreover, we 
obtain ((ST} thanks to the bound: || • < k~ 2 \\ ■ \\ 2 k y 

To obtain the estimate (|83|). we subtract ((71)) from l(77l) to obtain 



c fc ~~ C fc 



{ 



(/* I Ps) (f I P?) \ j f (4 I «?) (** I ft) 



INIo.i 



yiii2 

J s II 0.1 



llpfllg,, IMI§,i 



dcf 



We bound I\ by Lemma f5.ll < h 2a ° ||/ fe || 01 - As for I2, it is zero when [i = 0; 
otherwise we rewrite it as follows: 



h 



5T« 



iirfiiii 
+ 



+ ( 4 - zk I 9s " 9s 



z k I gj - g* 



)} 



llrtllS. IWI§. 



n 1 S 



— ./g + jf + jf + Jj^ . 



Then, recalling that ||g s || is constant, we estimate: 

• 141 £ 114 - *lo,i £ [^ 2ao fc2(Q0_1) + h2 } \\f%,i b y ©. 



I Jf I < ||jz* - z fc || Q _ 1 < [/i 2a ° fc Q "- 2 + A:" 1 ] 



lo,i 



by @g> and (EH). 
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• \4\ ^ h2ao \\ zk \\ s ~ h2ao k ~ 2 \\f k \\ 1 b y <E3> and <H3>- 

• | J 2 4 | < /i 2q o ||z fc || 1 < h 2a ° r max ||2 fe || _ x < h 2a ° k- 2 ||/ fe || 0jl by @B and JTSJ). 
Summarising, we obtain 



< (h 2ao + /i 2a ° + h 2a ° k 2ao + /i 1+Q ° fc + h 2 k 2 ) 



f 



0,1 



The estimate (|83j) then follows by remarking that the first, second and fourth terms in 
the bracket are negligible with respect to the third. o 

Now, we observe that the formula (|79|) for computing u%, as well as the SCM 
reconstruction formula for u k 



«£ = «* + 4 <Ps=^h + 4 {fa + s h <p P ). 



h r,% 



(84) 



are formally identical to their prismatic counterparts (cf. the SCM algorithm of |12l 

§5.3]); and the "building blocks" 4 an d a l so satisfy estimates similar to those of 

i 

the prismatic case. Indeed, under the assumption k < C* h 2 ~ a o , both terms within 
the bracket in (|83[) are negligible with respect to hk. Hence the following two results, 
whose proofs closely parallel that of Lemma 5.4 and Theorem 5.1 in ^2] 5 with the same 
kind of adaptations (use of weighted norms, P/j and Lemma l4.2j) as usual. 



Lemma 5.4 The following error estimate holds: 

2 



4 



(*) 



< k { h (1 + k h 



2 u2\ 



+ 



0,1 



ck-4 



Theorem 5.1 Let u k be the solution to the equation and u k be its finite element 

approximation given in \8J\) . Then the following error estimate holds: 



u k -u k h 



(k) 



< k 2 h 



f 



0,1 



(85) 



5.4 Approximation of the singular functions for the modes 

|fc| = 0, 1. 

The FE approximation of these functions has been exposed in [SJ §§4.1 and 4.2]. (In 
that work, the Laplacians Ao and Ai are respectively called A and A'). We keep 
this method, with the following modification. The dual singular functions associated 
to the reentrant edge undergo the same inconvenience as p 2 , namely, the Laplacian of 
the principal parts as defined in [2] do not belong to L\(uj). Hence we must enrich 
them, just as we did for p 2 p , in order to preserve the convergence rate. Calculating like 
in H5.1( we obtain the following decompositions: 



1 1 . . — 1 1 def _ a . / ,% f 

P s =Pp+P , P p = P sm(a0) - 
F 1 a 



1 3 P jj 
1 cos 4> 



P. 



0,e 



n e def _ a . / 

Pp = p sin(a(j 



2a 
P_ 

2rt 



COS i 



?€V\(u); (86) 
pP' e £ H\(u). (87) 
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The Laplacians of the principal parts are: 



il 4£. f 



„qO,e def . n,e 



1 



9 15 

- sin(a0) + — a cos sin(a0 + <j)') 

13 

-— sin(a0) + — a cose/)' sin(a0 + </> ) 



G L?( W ) 

e L\{u). 



Then we proceed like in H5.ll to obtain: 





pI-pI* 


1,1 
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0,e 0,e;h 
s Ps 












1,1 







Ps 



0.1 



< h 



2a i 



0.1 



(89) 

(90) 
(91) 



However, for the primal edge singular functions, the method of 5 yields the desired 
convergence rate. We just recall the decompositions: 



as well as the Laplacians of the principal parts: 

sin [(a - 1) 



<P k = (£) P° sin(a<5 



,fcdef A fc (fc + l)(fc + 2) 



a p 



2a fc r 1 k 

The line of proof already exposed in H5. 21 then easily leads to the error estimates: 



Wl - ^lli i ~ h and M 



,0,e 
s 



09°' e;/l | </l 
r^s ll 1 ~ 



(92) 



(93) 



(94) 



Now, as far as the conical point singularities are concerned, the method appears 
very similar to that of |12l §§5.1 and 5.2] since the principal parts p°^ c and <^P' C have a 
vanishing Laplacian Aq. So, mutatis mutandis, we get the error estimates: 



vY - Pi c '" 





















<h. 



i.i 



(95) 



Remark 5.1 Thanks to the asympotic expansions Eq. 8.7.1] of the Leg endre func- 
tion, it is possible to compute the function P u (cos (f>) with an arbitrary precision. Thus, 
one can compute once and for all the singularity exponent v and the integral in \^Ujj 
with an accuracy equal to the machine precision. All this guarantees that the errors 
due to the approximation of the conical singular functions will be negligible before the 
FE discretisation error. 



5.5 Approximation of u k and c/-, for \k\ < 1. 

As the representation formulae ()33|) and Q46|) for the singularity coefficients of these 
modes are rather standard, one can use the simple discrete versions: 

/.-til l:h\ ( rQ i 0,j;h 

f IPs' , / U I P. 



C ±l = il l;fei|2 > C 0. = || 0,y,H ]]2 J ■ (96) 

II II 0,1 WP S 1 1 0,1 
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Similarly, we will approximate the regular parts u k , \k\ < 1, by ul, u, £ V^, and 



iiP h G V; ft such that 



k = ±1 : ai u fc ) + 4 ai u & ) = (/* | , Vv A 6 , 

= : a (4, Wfc) + cj e «0 (<P° s ' e ' h , v h ) + c^ c a (</# c;ft , ifc) = (/° | « ft ) , to A € V ft 
Of course, we have the SCM reconstruction formulae: 



k = ±l: u k h = u k h + 4 ifl' h =u k h + 4 + sW P ) ; 



r,0 



i h I ~0,e . r-u.e u,e\ , h ( ~u,c . ru.c i 
u h + C , e (ifh + <V ^ ) + c , c + S h ' 

The results of £15. 41 then allow to conclude that: 



s; >e 0,. 



Ji f.~0,c 



°> c ,„0,c\ 
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c fe 



< h 



l ufc -^lln 



<h 



0,1 
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0.1 



c o,j - c o,j 



<h\\f 



1 0,1 



\u°-u° 



°l < ft II f°ll 
^li,i ~ Ik 1 1 o,i ' 



|„° _ ,.0| < u II f0\\ 
\ U U h\l,l ~ A ||/ 1 1 0,1 



(99) 
(100) 



(101) 
(102) 



6 Fourier Singular Complement Method 

Let u be the solution to the 3D problem (|2"|). and its Fourier coefficients. From 
the previous Sections, we know that u k (r,z) solves the 2D problem Q, the weak 
formulation of the elliptic problem (JJHSJ). And, according to the mode k, one can 
decompose u k as ([T7|l. (J32J) or (gSJ). 

The result of Heinrich |15l Thm 5.2] can be straightforwardly extended to our 
domain with a sharp vertex. 

o 

Theorem 6.1 Let f G h 2 (Q), and u € H 1 ^) be the solution to HJ). Then: 

u(r, 9, z) = u(r, 6, z) + 7 (0) tf(r, z) + c c Q ip° s ' c (r, z), (103) 

O 

with: u G H 2 (Q) fl i? 1 ^), and 7 E ^ 2 (S 1 ) is given by the formula: 

7(«o = ^<« + £ E c * eifce + E c ^ 9 - 

fe=±l |fe|>2 

Like in the prismatic case (cf. |12j . Remark 6.1), the hypothesis / £ /i 2 (J7) is crucial: 
the lack of its satisfaction would prevent the convergence of 7 in a regular enough 
space, and hence that of the singular part of the solution in the natural space. 
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We define the Fourier-SCM (FSCM) solution to © as follows: 



u 



[N] 



N 



k=-N 



where is the SCM solution to @ algorithmically defined in £J5J The main result on 
this method is the following 

Theorem 6.2 Assume f E h 2 (£2). Then the following error estimate holds: 

I f 



U — U 



[N] 
h 



Proof: Using the definition of and © we have 



u — u 



[N] 
h 



u 



\k\<N 

Using (fT§)) . we estimate E2 as: 
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clef 



\k\>N 
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E 2 < N~ 2 E fc2 

|fc|>iV 



1.1 



7 |fc|>iV 



< JV 



-2 



0.1 



2 

L 2 (n) ■ 



As for £1, we cut it into three parts, corresponding to k = 0, \k\ = 1, and 2 < < iV, 
which we bound respectively by (|l()2j) and (jlOlj) and (j%5|) : 

2 



^ £ ^ 2 ll/t, 1 + ^ 2 (ll/ 1 !lo, 1 + Ilr 1 lli0+^ E fc4 

2<|fc|<iV 

5 2 / 



r 



0,1 



< ^ 



L 2 (H) + 



de 2 



L 2 (n) 



where we have used Lemma 12.21 to bound the sum. Hence the result. 



7 Conclusion 

In this paper, we have proven that the FSCM for the Poisson equation achieves the 
optimal convergence rate for Pi finite elements and a datum of L 2 -style regularity in 
the meridian directions. The same result also holds for the discretization of the Poisson 
problem with a homogeneous Neumann boundary condition, or with non-homogeneous 
boundary conditions, provided there exist sufficiently smooth liftings. 

This result closely parallels that of the companion paper ^5]. The specificities of 
the axisymmetric geometry (namely, that the 2D problems are set in weighted Sobolev 
spaces, which moreover vary for the low-order Fourier modes before stabilising, and 
involve differential operators with non-constant coefficients) only cause technical dif- 
ficulties. As far as the presence of conical vertices is concerned, its effect is no more 
than a finite-dimensional perturbation. Furthermore, it is no difficulty to consider the 
case of an axisymmetric domain Q with several reentrant edges (i.e. u with several 
off- axis reentrant corners) and/or several sharp vertices. 
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As already mentioned, this paper is the second part of a three-part article |12l 113] . 
In ^3], the FSCM is analysed from a numerical point of view (complexity, implemen- 
tation issues, numerical experiments, etc.), and it is compared to other methods — in 
the axisymmetric case, to anisotropic mesh refinement techniques. 

One can apply the same theoretical and numerical techniques to the fully axisym- 
metric heat or wave equations, with any L 2 -smooth (in space) right-hand side. For 
these PDEs, the singular functions p s and <p s do not depend on the time-step. 
Finally, the results, can also be viewed as the first effort towards the discretization 
of electromagnetic fields in axisymmetric domains, with continuous numerical approx- 
imations, the importance of which is well-known, cf. [2j. As a matter of fact, the 
SCM developed in [31 0J |S] for fully axisymmetric electromagnetic computations can 
be generalized to arbitrary data, with the help of the results obtained here. 

References 

[1] Abramowitz, M., Stegun, LA. (1965): Handbook of mathematical functions. Dover, New 
York. 

[2] Assous, F., Degond, P., Hcintzc, E., Raviart, P.-A., Segre, J. (1993): On a finite-element 
method for solving the three-dimensional Maxwell equations. J. Comput. Phys. 109, 
222-237. 

[3] Assous, F., Ciarlet, Jr, P., Labrunie, S. (2002): Theoretical tools to solve the axisymmetric 
Maxwell equations. Math. Meth. Appl. Sci. 25, 49-78. 

[4] Assous, F., Ciarlet, Jr, P., Labrunie, S. (2003): Solution of axisymmetric Maxwell equa- 
tions. Math. Meth. Appl. Sci. 26, 861-896. 

[5] Assous, F., Ciarlet, Jr, P., Labrunie, S., Segre, J. (2003): Numerical solution to the 
time-dependent Maxwell equations in axisymmetric singular domains: the Singular Com- 
plement Method. J. Comput. Phys. 191, 147-176. 

[6] Assous, F., Ciarlet, Jr, P., Sonnendriickcr, E. (1998): Resolution of the Maxwell equations 
in a domain with reentrant corners. Math. Model. Numer. Anal. 32, 359-389. 

[7] Assous, F., Ciarlet, Jr., P., Segre, J. (2000): Numerical solution to the time-dependent 
Maxwell equations in two-dimensional singular domains: the Singular Complement 
Method. J. Comput. Phys. 161, 218-249. 

[8] Belhachmi, Z., Bernardi, C, Deparis, S. (2003): Weighted Clement operator and appli- 
cation to the finite element discretization of the axisymmetric Stokes problem. Internal 
Report 03029, Laboratoire Jacques-Louis Lions, Universite Pierre et Marie Curie. 

[9] Bernardi, C, Dauge, M., Maday, Y. (1992): Polynomials in weighted Sobolev spaces: 
basics and trace liftings. Internal Report 92039, Laboratoire d'analyse numerique, Uni- 
versite Pierre et Marie Curie. 

[10] Bernardi, C, Dauge, M., Maday, Y. (1999): Spectral methods for axisymmetric domains. 
Series in Applied Mathematics, Gauthier-Villars, Paris and North Holland, Amsterdam. 

[11] Ciarlet, Jr, P., He, J. (2003): The singular complement method for 2D scalar problems. 
C. R. Acad. Sci. Paris, Ser. 1 336 353-358. 

[12] Ciarlet, Jr, P., Jung, B., Kaddouri, S., Labrunie, S., Zou, J. (2004): The Fourier Singular 
Complement Method for the Poisson problem. Part I: prismatic domains. To appear in 
Numer. Math. 



27 



[13] Ciarlet, Jr, P., Jung, B., Kaddouri, S., Labrunie, S., Zou, J. (2004): The Fourier Singular 
Complement Method for the Poisson problem. Part III: implementation issues. In prepa- 
ration. 

[14] Grisvard, P. (1992): Singularities in boundary value problems. RMA 22, Masson, Paris. 

[15] Heinrich, B. (1993): Singularity functions at axisymmetric edges and their representation 
by Fourier scries. Math. Meth. Appl. Sci. 16, 837-854. 

[16] Heinrich, B. (1996): The Fourier-finite element method for Poisson's equation in axisym- 
metric domains with edges. SI AM J. Numer. Anal. 33, 1885-1911. 

[17] Heinrich, B., Nicaise S., Weber, B. (2000): Elliptic interface problems in axisymmetric 
domains II: Convergence analysis of the Fourier-finite element method. Adv. Math. Sci. 
Appl. 10 571-600. 

[18] Mercier, B., Raugel, G. (1982): Resolution d'un probleme aux limites dans un ouvert ax- 
isymetrique par elements finis en r, z et series de Fourier en 9. RAIRO Analyse numerique 
4, 405-461. 

[19] Nazarov, S.A., Plamenevsky, B.A. (1994): Elliptic problems in domains with piecewise 
smooth boundaries. De Gruyter Exp. Math. 13. 

[20] Raugel, G. (1978): Resolution numerique de problemes elliptiques dans des domaines avec 
coins. PhD Thesis, Rennes University, France. 

[21] Zienkiewicz, O.C. (1977): The Finite Element Method. McGraw-Hill, London. 



28 



